#ifndef VEC_H_
#define VEC_H_

#include "../common.h"

#include "math_exception.h"
#include "aux_funcs.h"
#include <valarray>
#include <ostream>

namespace physicsworks {

namespace math {

using exception::EBadArg;
using exception::EUndefMathOpt;

template<class T> class Vec {
protected:
	std::valarray<T> arr; /**< the vector in low-level */

public:
	typedef T Type;

	Vec<T>(size_t sz) : arr(sz) {
		Assert<EBadArg>(!DMATH || sz > 0, "The size of a vector must be greater than 0");
	}
	Vec<T>(const Vec<T>& v) : arr(v.arr) {}
	
	size_t size() const { return arr.size(); }
	
	/**
	 * Returns the i-th element.
	 */
	T operator () (size_t i) const {
		Assert<EBadArg>(!DMATH || i >= 0, "The index must be non-negative");
		Assert<EBadArg>(!DMATH || i < arr.size(), "The index must be smaller than the size");
		
		return arr[i];
	}
	
	/**
	 * Returns the i-th element.
	 */
	T& operator () (size_t i) {
		Assert<EBadArg>(!DMATH || i >= 0, "The index must be non-negative");
		Assert<EBadArg>(!DMATH || i < arr.size(), "The index must be smaller than the size");
		
		return arr[i];
	}
	
	/**
	 * Assigns the other vector.
	 */
	Vec<T>& operator = (const Vec<T>& v) {
		Assert<EUndefMathOpt>(!DMATH || v.size() == size(), "Both vectors must have the same size");
		
		arr = v.arr;
		return *this;
	}

	/**
	 * Sums by other vector.
	 */
	const Vec<T> operator + (const Vec<T>& v) const {
		Vec<T> r = *this;
		return r += v;
	}
	
	/**
	 * Subtracts by other vector.
	 */
	const Vec<T> operator - (const Vec<T>& v) const {
		Vec<T> r = *this;
		return r -= v;
	}

	/**
	 * Multiplies by other vector.
	 */
	const Vec<T> operator * (const Vec<T>& v) const {
		Vec<T> r = *this;
		return r *= v;
	}
	
	/**
	 * Divides by other vector.
	 */
	const Vec<T> operator / (const Vec<T>& v) const {
		Vec<T> r = *this;
		return r /= v;
	}

	/**
	 * Sums by a scalar.
	 * All vector' scalars are summed by the provided scalar.
	 */
	const Vec<T> operator + (T v) const {
		Vec<T> r = *this;
		return r += v;
	}
	
	/**
	 * Returns a new vector by applying minus on all elements.
	 */
	const Vec<T> operator - () const {
		Vec<T> r(arr.size());
		r.arr = -arr;
		return r;
	}
	
	/**
	 * Subtracts by a scalar.
	 * All vector' scalars are subtracted by the provided scalar.
	 */
	const Vec<T> operator - (T v) const {
		Vec<T> r = *this;
		return r -= v;
	}

	/**
	 * Multiplies by a scalar.
	 * All vector' scalars are multiplied by the provided scalar.
	 */
	const Vec<T> operator * (T v) const {
		Vec<T> r = *this;
		return r *= v;
	}
	
	/**
	 * Divides by a scalar.
	 * All vector' scalars are divided by the provided scalar.
	 */
	const Vec<T> operator / (T v) const {
		Vec<T> r = *this;
		return r /= v;
	}

	/**
	 * Sums by other vector and assigns.
	 */
	Vec<T>& operator += (const Vec<T>& v) {
		Assert<EUndefMathOpt>(!DMATH || v.size() == size(), "Both vectors must have the same size");
		
		arr += v.arr;
		return *this;
	}
	
	/**
	 * Subtracts by other vector and assigns.
	 */
	Vec<T>& operator -= (const Vec<T>& v) {
		Assert<EUndefMathOpt>(!DMATH || v.size() == size(), "Both vectors must have the same size");
		
		arr -= v.arr;
		return *this;
	}

	/**
	 * Multiplies by other vector and assigns.
	 */
	Vec<T>& operator *= (const Vec<T>& v) {
		Assert<EUndefMathOpt>(!DMATH || v.size() == size(), "Both vectors must have the same size");
		
		arr *= v.arr;
		return *this;
	}

	/**
	 * Divides by other vector and assigns.
	 */
	Vec<T>& operator /= (const Vec<T>& v) {
		Assert<EUndefMathOpt>(!DMATH || v.size() == size(), "Both vectors must have the same size");
		
		arr /= v.arr;
		return *this;
	}

	/**
	 * Sums by a scalar and assigns.
	 * All vector' scalars are summed by the provided scalar.
	 */
	Vec<T>& operator += (T v) {
		arr += v;
		return *this;
	}
	
	/**
	 * Subtracts by a scalar and assigns.
	 * All vector' scalars are subtracted by the provided scalar.
	 */
	Vec<T>& operator -= (T v) {
		arr -= v;
		return *this;
	}

	/**
	 * Multiplies by a scalar and assigns.
	 * All vector' scalars are multiplied by the provided scalar.
	 */
	Vec<T>& operator *= (T v) {
		arr *= v;
		return *this;
	}

	/**
	 * Divides by a scalar and assigns.
	 * All vector' scalars are divided by the provided scalar.
	 */
	Vec<T>& operator /= (T v) {
		arr /= v;
		return *this;
	}
	
	bool operator == (const Vec<T>& v) const {
		for (size_t i = 0; i < arr.size(); ++i)
			if ( !equal(arr[i], v(i)) )
				return false;
		return true;
	}
	
	bool operator < (const Vec<T>& v) const {
		for (size_t i = 0; i < arr.size(); ++i) {
			if ( less(arr[i], v(i)) )
				return true;
			else if ( greater(arr[i], v(i)) )
				return false;
		}
		return false;
	}
	
	/**
	 * Returns an array of reals representing the vector.
	 */
	real* toReal() const {
		real* r = new real[arr.size()];
		for (size_t i = 0; i < arr.size(); ++i) r[i] = arr[i];
		return r;
	}
	
	void zero() {
		for (size_t i = 0; i < arr.size(); ++i) arr[i] = 0;
	}
};

}

/**
 * Print the VecBasis<T> formatted.
 */
template<class T> inline std::ostream& operator << (std::ostream& os, const math::Vec<T>& v)
{
	os << "Vec(";
	size_t i = 0;
	goto entry;

	do {
		os << ", ";
entry:
		os << v(i++);
	} while(i < v.size());
	
	os << ")";
	return os;
}

/**
 * Scalar as left-side version from operator*(T v).
 */
template<class T> inline math::Vec<T> operator * (T l, const math::Vec<T>& v)
{
	return v * l;
}

/**
 * Scalar as left-side version from operator+(T v).
 */
template<class T> inline math::Vec<T> operator + (T l, const math::Vec<T>& v)
{
	return v + l;
}
	
}

#endif /*VEC_H_*/
